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Abstract 


This  paper  presents  an  algorithm,  the  secant/finite  difference  algorithm,  for  solving  sparse 
nonlinear  systems  of  equations.  This  algorithm  is  a  combination  of  a  finite  difference  method  and 
a  secant  method.  A  5-superlinear  convergence  result  and  an  r-convergence  rate  estimate  show 
that  this  algorithm  has  good  local  convergence  properties.  The  numerical  results  indicate  that  this 
algorithm  is  probably  more  efficient  than  some  currently  used  algorithms. 


1.  Introduction. 


Cionsider  the  nonlinear  system  of  equations 

F(i)  =  0,  (11) 

where  F:f?*  —►/?“  is  continuously  differentiable  on  an  open  convex  set  C  F",  and  the  Jacobian 

matrix  F  '(z)  is  sparse.  To  solve  the  system,  we  consider  the  iteration 

X  =x-B-‘F(*).  (1-2) 

where  x  is  the  current  iterate,  x  is  the  new  iterate  and  B  is  an  approximation  to  the  Jacobian 
F  “’(z)  with  the  same  sparsity  as  the  Jacobian.  After  we  finish  this  step,  we  have  the  information. 
X,  z,  and  B.  Our  purpose  is  to  get  a  B,  a  good  approximation  to  F'(z),  by  using  as  little  effort 
as  possible. 

Curtis,  Powell  and  Reid  [3]  gave  an  efficient  algorithm  for  sparse  problems  called  the  CPR 
Algorithm,  which  is  a  finite  difference  algorithm,  but  which  can  take  the  advantage  of  the  sparsity 
to  make  the  number  of  the  function  evaluations  small.  To  describe  the  CPR  Algorithm,  Coleman 
and  More  [2]  gave  some  definitions  concerning  a  partition  of  the  columns  of  a  matrix  B. 

Definition  1.1.  A  partition  of  the  columns  of  B  is  a  division  of  the  columns  into  groups 
Cl,  Cj  •  •  •  Cf  such  that  each  column  belongs  to  one  and  only  one  group. 

Definition  l.S.  A  partition  of  the  columns  such  that  columns  in  a  given  group  do  not  have  a 
nonzero  element  in  the  same  row  position  is  consistent  with  the  direct  determination  of  B. 

In  order  to  have  a  good  partition  of  the  columns  of  B,  Coleman  and  More'  [2]  associated  the 
partition  problem  with  a  certain  graph  coloring  problem  and  gave  some  partitioning  algorithms 
which  can  make  p  optimal  or  nearly  optimal.  For  convenience,  we  call  the  CPR  algorithm  based 
on  Coleman  and  More'’s  algorithms  the  CPR-CM  algorithm. 

The  CPR  Algorithm  can  be  formulated  as  follows:  For  a  given  consistent  partition  of  the 
columns  of  B,  obtain  vectors  d,,  d^,  ■  ■  ■  d^  such  that  B  is  determined  uniquely  by  the  equations 

Bdi  =  F(x  +  d.)  -  F(x)  =  Vi,  i  =  1,  P  ■ 


Notice  that  at  each  iterative  step  we  need  only  to  compute  p  +  1  function  values  rather  than  the 
n  +  1  values  required  by  a  straightforward  column-by-column  finite-difference  algorithm. 

As  an  example  we  consider  the  matrix  with  a  tridiagonal  structure: 

X  X  0  0  0  o' 

X  X  X  0  0  0 

0  X  X  X  0  0 

OOXXXO- 
0  0  0  X  X  X 

.0  0  0  0  X  X. 

A  consistent  partition  of  the  columns  of  the  matrix  is  Cj  =  {l,  4},  C2  =  {2,  5},  and  Cg  =  {3,  6}. 
If  we  take 

rfi  =  (A  ,0,0,A  ,0,0)^, 
da  =  (  0  ,  A  ,  0 , 0  ,  A  ,  0  )^, 
d3  =  (0,0,  A  ,0,0,  A 

then  B  is  determined  uniquely,  and  the  number  of  function  evaluations  required  at  each  iteration 
is  4. 

In  this  paper,  we  propose  an  algorithm  called  the  secant/finite  difference  (SFD)  algorithm 
for  solving  sparse  nonlinear  systems  of  equations.  This  algorithm  is  also  based  on  a  consistent 
partition  of  the  columns  of  the  Jacobian.  However,  it  uses  the  information  we  already  have  at 
every  iterative  step  more  efficiently  than  the  CPR  algorithm.  The  secant  equation  is  also  satisfied 
by  the  SFD  algorithm.  Therefore,  this  algorithm  can  be  seen  as  a  combination  of  the  CPR-CM 
algorithm  and  a  secant  algorithm.  The  SFD  algorithm  reduces  the  number  of  function  evaluations 
required  by  the  CPR-CM  algorithm  by  one,  and  it  has  good  local  convergence  properties.  Our 
numerical  results  show  that  the  SFD  algorithm  is  probably  more  efficient  than  the  CPR-CM  algo¬ 
rithm.  The  SFD  algorithm  and  some  of  its  properties  are  given  in  Section  2.  A  Kantorovich-type 
analysis  for  the  SFD  algorithm  is  given  in  Section  3.  A  g-superlinear  convergence  result  and  an 
r-convergence  order  estimate  of  the  SFD  algorithm  are  given  in  Section  4.  Some  numerical  results 
are  given  in  Section  5. 
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In  thix  paper,  ll.Hf  indicates  the  Frobenius  norm  of  a  matrix,  and  ||.||  indicates  the  /jr 
vector  norm.  We  use  ‘  \  ’  to  denote  the  subtraction  of  two  sets;  that  is, 

A  \  B  =  {v  :  V  E  A  and  v  ^B}  . 

For  a  sparse  matrix  B,  we  use  M  to  denote  the  set  of  pairs  of  indices  (i,  j),  where  6,'y  is  a  struc¬ 
turally  noniero  element  of  B,  i.t. 

M  =  j)  :  bi,  0}  . 


2.  The  SFD  Algorithm  and  its  Properties. 

Given  a  consistent  partition  of  the  columns  of  the  Jacobian,  which  divides  the  set  {1, ...,  n} 
into  p  subsets  ei,...,Cp  (for  convenience,  c,-,  i=l,2,...,p,  indicates  both  the  sets  of  the  columns 
and  the  sets  of  the  indices  of  these  columns),  also  given  x,  x  €  J?",  let 

d,-  =  ^  ]  SyCy,  I  1,  *  P  f 
see, 

and 


Qi  —  ^  ^  dy ,  I  —  1, ...,  p  ,  ^0  * 

y-1 

If  «y  /  =  1,  ...n,  then  B  is  determined  uniquely  by  the  equations 
Bdi  =  .B(x-(x-f,))  =  F{x)-F{x-gi)  =  pi  , 

Bdi  =  .B(x-ff,-.i-(x-ji,))  =  F(x-y,_i)-F(x-f,)  =  y,-  ,  (2.1) 

Bd^  =  B{x-g,_i-x)  =  F{x-gp_i)-F{x)  =  y,  . 

Let  B  =  [6jb,].  By  (2.1),  if  (/,  m)  E  M,  then 


where  tnECf,  i—l,  2,  ...,  p. 

Notice  that  by  (2.1),  to  get  B,  we  need  only  to  compute  p  function  values  since  we  already 
have  the  value  F(x)  at  the  current  step.  The  number  of  function  evaluations  at  each  iteration  is 


S 


one  less  than  the  CPR-CM  algorithm.  For  example  (1.3),  we  take 

rfi  =  (  ,  0 , 0  ,  <4 , 0 , 0  )^, 

rfj  =  (  0  ,  <2  >  ®  ®  »  *6  >  0  )  > 

~  (  0  >  0  >  •$  >  ®  ®  >  *8  )  > 

and 

,0,0,84,0,0)^, 

92  (  8i  ,  82 , 0 , 84  ,  85 , 0  ) 

Then,  we  need  only  to  compute  the  values  of  F{x),  F(x— ji),  and  F{x—g^,  and  the  number  of 
function  evaluations  is  3  per  iteration  instead  of  4  required  by  the  CPR-CM  algorithm. 

Let 


Ji=^  i  F'{x-gi  +  tdi)it,  •=l,...,p 


(2.3) 


Then 

Jiii^Vi,  1  =  1,  ...,p.  (2-4) 

Let  /,•  =  [Jim]-  Since  J,-  has  the  same  sparsity  as  the  Jacobian,  by  (2.4),  we  have  that  if 


(/ ,  m )  €  M,  then 


T 

Cl  Vi 


•>lm 


where  mEe,-,  «=1,  2,  ...,  p.  Comparing  (2.5)  with  (2.2),  we  have 

J,*  Cj  —  Bcj  , 

where  j  €  c,-  ,  »  =  1, P  •  Therefore,  B  can  be  written  as 

•  =  1  i€c. 


(2.5) 


-(2.6) 


(2.7) 


To  study  the  properties  of  the  SFD  algorithm,  we  assume  that  F  (x )  satisfies  the  Lipschitz 
condition,  i.e.  that  there  exist  a,-  >  0,  i  =  1, ...,  n  such  that 

||(F'(x)-F'(y))e,||  <  a.lk-yll,  f  =  l,2,...,n  ,  x,yED.  (2.8) 


Let  a  =  ( 5^ 

f  =  l 
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|lF'(x)-F'(»)l|r  <alU-|r|l.  x,ir€D. 

Now  we  have  the  following  estimate  for  S : 

Theorem  2.1.  Suppose  F"{x)  satisfies  Lipschitz  condition  (2.8),  and  B  i 
xeD,  x-giED,  i  =  l,...,p,  and  «,  7^  0,  i=l,  2,  n,  then 

||F'(x)-B||f  <  a|l*-xl|  . 

Proof.  By  (2.6)  and  (2.7), 

||F'(x)-Bl|?=  E  |1(F'(x)-B)c,||2 
•—1 

=  £  E  \\{F'ix)-B)e,-r 

•  -=1  iee, 

=  £  E  \\{F’{x)-J,)e,\\^  . 

•  ■=1  i€c. 

Using  (2.3)  and  Lipschitz  condition  (2.8),  we  obtain 

E  \\{F'ix)-Ji)e,.\\^ 

1 

=  E  WiF'ix)-!  F'{x-gi  +  t{gi-9,_i))dt)e,  \\^ 

jet,  0 

<  E  («i/ 

iec,  0 

1  1 

<  E  a/(/(l-OIU.II‘i<  +  /lU.-ill‘‘^‘f 

}€c,  0  0 

i€c,  ^  ye«, 

It  follows  from  (2.10)  and  (2.11)  that 

\\F'{x)-B\\j<  ||.|r  £  E  of=a^\\s\\ 

•==i  yec, 


determined  by  (2.1).  If 

(2.9) 

(2.10) 


(2.11) 

(2.12) 
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Then  (2.9)  follows  from  (2.12). 

In  (2.1),  to  determine  S  uniquely,  we  assume  that  Sy  j  =  1, ...,  n .  However,  sometimes 
it  may  happen  that  «,  =0  for  some  1  <  »  <  n.  If  this  happens,  then  the  ith  column  of  B  can 
not  be  determined  uniquely  by  (2.1).  In  this  case,  let 

=  {»■  G  {1,  2,  ...,  n):  «,  0)  , 

and  let 

n2={l,  2,  ...,  n)\n,  . 

Now  we  deal  with  the  general  case  in  such  a  way  that  if  j  ECly,  then  the  yth  column  of  B  is 
determined  uniquely  by  (2.1).  If  j  €^2,  then  we  let  the  jth  column  of  B  be  equal  to  the  jth 
column  of  B.  In  practice,  if  l«y|  is  too  close  to  zero  the  cancellation  errors  will  become 
significant.  Therefore,  there  should  be  a  lower  bound  0  for  j  «y  | .  Now  the  SFD  algorithm  with  a 
global  strategy  can  be  stated  as  follows; 


Algorithm  2.S.  Given  a  consistent  partition  of  the  columns  of  the  Jacobian,  which  divides  the  set 
{1,  2,  ...,  n)  into  p  subsets  Cj,  C2,  ...,  Cp,  and  given  ar  ‘  ,  x®  €  B"  such  that 
s  x°  -  0  >  2,  •  •  •  ,  n,  at  each  step  k  >  0: 

(1).  Set 

=  E  E  .  '=1-  2. P-i ,  go-^  =  o  • 
y=i 


where  =  x*  -  i*  ^ 

(2) .  Ckimpute  B(x*  -  ji*"*),  i=0,  1,  ...,  p-1,  and  set 

=  F(x*  -  )  -  B(x*  -  gr)  .  •••.  P 

where  F(x*  -  H*'*)  =  B(x*"')  . 

(3) .  If  (/,  m)  €  M  and  |  «*"*  |  >  9,  then  set 

,  r,,*-i 

A*  —  ‘ 

otherwise  set 


(2.13) 
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ik  _  JL*-1 

Olm  —  ®l«  . 

where  m€e,-  ,*=1,  2,  p  . 

(4) .  Solve  ==  -F(z*). 

(5) .  Choose  by  x*+*  =  x*  +  «*,  or  by  a  global  strategy. 

(6) .  Check  for  convergence. 


The  SFD  algorithm  is  also  an  update  algorithm,  and  the  update  can  be  written  as 


yefij  ••=1  yec.nOi 

The  following  result  shows  that  the  SFD  algorithm  is  a  secant  algorithm. 
Lemma  £.S.  B  satisfies  the  secant  equations 

Bd{  =  y{,  i  =  l,...,p  , 


(2.14) 


(2.15) 


and  (2.15)  implies 

Ba  =  F(I)-F(x)  =  !/  . 

The  proof  of  Lemma  2.3  is  straightforward. 

Suppose  that  we  have  finished  the  jfcth  step  of  the  iteration.  Then  the  information  we  have 
is  X*,  F(x‘),  Bt,  and  x*+*.  Let 

=  E  >  «■  =  1-  •  •  •  P  > 

yec, 


?*=E<^L  1  =  1,...,?,  Jo=0, 

y=i 


and 


=  /  F'(x*+‘-p,^+td,^)dt ,  •  =  l,...,p.  (216) 

0 

Then  by  (2.6), 

=  Bt+iCy  ,  (217) 

where  j  6c,-  ,  »=1,  2,  ...,  p. 

Theorem  2.4.  Assume  that  F'  satisfies  Lipschitz  condition  (2.8).  Let  {x'}*±o  {^y}y=o  be  gen¬ 
erated  by  Algorithm  2.2.  Suppose  that  i  =  0, 1, 2,  ...,p }y=o  C  D.  If  «,*  =  0  appears 
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consecutively  in  at  most  m  steps  for  any  specific  1  <  i  <  n ,  then  for  t  >  m , 


<  o  E  (2.18) 

y»*-m 

Proof.  By  the  hypothesis  of  the  theorem,  given  i,  for  any  1  ^  w,  there  exists  at  least  one 
integer  0  <  y  <  m  such  that  Let  j{k,i)  be  the  smallest  one  of  these  integers.  Then, 

Let  i€C(.  Then, 

by  (2.17).  Therefore, 

||(F'(x*+‘)-B*,i)«.ll 
=  ||(F'(:r*+>)-5*.y(*,)+i)«.|| 

<  ll(F'(x‘+')-F'(x*-^(*-‘»))c,||  +  ||(F'(x*--’(‘-*>0-5w(M)+i)«.  II 

=  ||(F'(x*+‘)-F'(x*-''(*'’)+‘))e.|l  +  ||(F'(x*-''(*''>+‘)-Jf-''(*'‘‘)^^^ 

<  a,- 1 1  x‘+‘  -  x*-^(*'‘)+*  1 1  +  a,- 1 1  x*-^(*'')+‘  -  x*-^'(*-’)  1 1 

<a,  E  ||x'+‘-x‘|| 

l^k-S(k,i) 

<a,  E 

tsssk-m 

Hence, 

||F'(x*+‘)-5*+iII?=  E  ll(F'(x‘+>)-F*+,)‘.ir 

(2.19) 

<(E 

Then,  (2.18)  follows  from  (2.19). 

3.  A  Kantorovich-Type  Analysis. 

The  following  estimate  for  the  SFD  algorithm  is  sharper  than  that  for  Broyden’s  algorithm 
given  by  Dennis  [4]. 
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Theorem  S.l.  Assume  that  F'  satisfies  Lipschitz  condition  (2.8)  and  that  {**}  and  {B*}  are  gen¬ 
erated  by  Algorithm  2.2  with  ||i~‘-z°||  <h  .  If  «  =0, 1, ..., p}*_o C  Z),  then 


<  oE  .  (3.1) 

y-o 

Proof.  Inequality  (3.1)  can  be  obtained  immediately  by  setting  m  =  i:-f  1  in  (2.18). 

Theorem  S.i.  Let  F'  satisfy  Lipschitz  condition  (2.8).  Suppose  that  x~^,  x^ED,  and  that  Bq, 
generated  by  i'*  and  z°,  is  a  nonsingular  n  X  n  matrix  such  that 


\x-^-x°\\<6,  |1Bo-M1f<A  l|Bo-‘i^(*“)ll  <»;, 


h  = 


<1, 


and 


{1-Za^6f  -  6  ’ 


o^<j  . 


(3.2) 


If  S{x°,2t*)  C  D,  where 


(3.3) 


then  {**},  generated  by  Algorithm  2.2  without  any  global  strategy,  converges  to  x*,  which  is  the 
unique  root  of  F(x)  in  S(x°,7)nZ),  where 


—  1  -  a/3S 


1+1 


(l-ajSSfJ 


Proof.  Consider  the  scalar  iteration 

=  ^f{tk),  to  =  0,  *=0,1,2,  •  •  •, 

where 


It  is  easy  to  show  that  the  sequence  { t*  }  satisfies  the  difference  equation 


From  this  equation,  it  follows  that  {4}  is  &  monotonically  increasing  sequence  and 


lim  tj,  =  t'  , 

A  — *  00 

where  t*  is  the  smallest  root  of  /(/)• 

Now,  by  induction  we  will  prove  the  following  estimate: 

||**+‘-x*||  <  <*+!-<*  .  t=0,l,2,  •  •  •  . 

For  ib  =  0,  we  have 

<  »?  =  ‘i-<o  • 

Suppose  that  (3.4)  holds  for  k  =  0, 1, 2,  ...m  —  1.  Then 

Therefore,  {x*}  C  5(x°,  <*),  and 

{x*-j^  ,=1,  •  •  •  p}  C  5(x°,2t*),  Jt=0,l,...,m. 

Using  Theorem  3.1,  we  have 

||B«-Bo||  <  I \B„-F'(z'”)\ U+  1  |F'(z-”)-F'(x“)|  If  +  1  |F'(x‘>)-5ol|f 

<2a53  I  I x'‘'''-x*  II  +  2a6  <  2a<*  +  2a{  <  . 

•  =o  ^ 

Then  by  Banach  lemma, 

Hence, 

I  lx’"+‘-x’"  II 

<  I  I  B-MI  F  I  m*’"  )-•?’(*’"■*) -*"-*)  I  I 

<  3^[-2.||x”-x’"->||  +«E*  I  |x*+‘-x*||  +a«]|  |x„-x„_i|| 

2  ,.=0 

—  K^m  ~  ^m-l)  +  ®  ^m-l  +  ~  ^m-l)  ~  ^m+l~^tn  • 

This  completes  the  induction.  By  (3.4),  it  is  easy  to  show  that  there  is  an  x*  €/)  such  that 

lim  X*  =  X*  . 

k—^oo 


(3.4) 
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The  uniqueness  of  x*  in  S(x°,  t  )ni?  can  be  obtained  from  Ortega  and  Rheinboldt’s  Theorem 
12.6.4  [8,  p.425]  by  setting  A(z)  ^  So- 

4.  Local  Convergence  Properties. 

To  study  the  local  convergence  of  the  SFD  algorithm,  we  assume  that  F:J?  C  S"  — has 
the  following  property: 

There  ie  an  x*  ED  ,  aueh  that  F(x*)  =  0  and  F'{x‘)  it  noneingular.  (4.1) 

Theorem  4.1.  Assume  that  F:D  C  S*  — ►S"  satisfies  (4.1)  and  that  F'  satisfies  Lipschitz  condi¬ 
tion  (2.8).  Let  {x*}  be  generated  by  Algorithm  2.2.  without  any  global  strategy.  Then  there  exist 
€,  5  >  0  such  that  if  x“*,  x°£D  satisfy 

l|x°-x*||  <  €,  llx-‘-x°||  <  6, 

then  {x*}  is  well  defined  and  converges  g-superlinearly  to  x*. 

Proof.  Notice  that  we  can  choose  c  small  enough  so  that  [  |  Bo^F{x°)  \  \  is  also  small  such  that 

h  <  ■^  and  that  5(x°, 2t*)CS,  where  h  and  t*  are  defined  by  (3.2)  and  (3.3)  respectively. 
6 

When  6  is  small  enough,  we  also  have 

where  is  defined  in  Theorem  3.2.  Therefore,  by  Theorem  3.2, 

{x*+‘-|- f  =  1,2,  ...,p}  C  ,  t=0,  1,  •  •  •  . 

Thus,  from  (2.14)  and  the  proof  of  Theorem  2.1,  we  have 

\\F'{x-)-B\rF 

=  EII(^'(=^')-^KII“+  E  ||(F'(x*)-b)c,||2 

I'en,  leOj 

=  £  E  ||(F'(x')-J,)c,.||2+  x;  ||(F'(x')-B)e,||“ 

1=1  yee.rmj  iea^ 


11 


=  £  E  11/  [F'{x*)-F\x-9i  +  t(g:-ti.,))]t,dt\\^ 

1-1  yee.nOj  0 

+  E  |l(F'(x*)-B)e,|r 

•eng 


<  a^\\x--J\\  +  \\x-x\\f+  ||F'(x*)-B||^ 
<a\2\\x--x\\  +  ||x*-x|ir+  l|^’'(x*)-B||^ 


Therefore, 

\\F'{x’)-B\\f  <  ||F'(x*)-B|lr  +  ZacT(x,x), 

where 


<7(x,x)  =  max{||x-x*l|,  ||x-x'l|}. 


Notice  that  by  Theorem  2.1  and  Lipschitz  condition  (2.8), 


<  ||F'(x*)-F'(x°)||f  +  ||F'(x°)-Bo|U 

<  a|lx'  -  x°||  +  a||x°-  x-‘|| 

<  a(€  +  5) 

Thus,  by  Dennis  and  More’s  [5]  Theorem  5.1,  we  know  that  {x*}  converges  at  least  g-linearly  to 


According  to  Dennis  and  Morels  [5]  Theorem  3.1,  to  get  g-superlinear  convergence,  we  need 
only  to  prove  that 


||(F*-F'(x*))(x*+»-x*)||  ^ 

||x*+‘-x*|l 


(4.2) 


If  for  all  1  <  « <  n,  «*  =  0  appears  consecutively  in  at  most  m  steps,  then  by  Theorem  2.4  it 
is  easy  to  show  that 


lim  1|F*-F'(x*)|l^=0. 

k-*oo 


(4.3) 
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Thus,  (4.2)  follows  immediately  from  (4.3). 

Otherwise,  let 

ylj  =  {i  e  {1,  2,  ...,  n}  :  For  any  k>0  ,  there  exute  at 
leaat  one  integer  m>k  such  that  s,"  7^  0  }  , 

and  let  j4i  =  {1, ...,  n  }  \  A2.  Then 

fi*-F'(x*)=  E  (B*-F'(x‘))e.«.^+  E  {Bi-F'ix‘))eie^  . 

•  €>4j  <£>42 

From  the  definition  of  Ai,  there  exists  a  large  integer  ko  such  that  «,*  =  0  for  all  t'GAi  and  k>ko. 
Therefore, 

E  {Bi,-F'{x‘))eie,^{x^*^-x^)  =  0  ,  (4  4^ 

«€/i 

for  k  >  ko-  Now  we  show  that 

limll  E  (B*-F'(x'))e.e,niF=0.  .45) 

*-*“  <€^2 

In  the  first  part  of  the  proof,  we  proved  that  lim  1|  x*-x'  |1  =  0.  This  implies  that  given  e  >  0, 

*-•00 

there  exists  an  integer  K  such  that 

By  the  definition  of  A^,  there  exists  an  integer  Ki,  which  depends  on  K,  such  that  for  every 
•  €  A2,  there  exists  at  least  one  integer  0  <  j  <  Ki  such  that  ^  0.  Let  K  =  K  +  Ki-  For 
k  >  K  and  1GA2,  define 

j{k,i)  =  min{  j  :  0). 

Then  k  -  i{i,k)  >  K.  Let  i€c,  ,  1  <  /  <  P,  we  have  that 

B,e,  =  . 

Thus,  by  Lipschitz  condition  (2.8), 
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1 

=  II  /{F  -  F'{x 

0 

<  (/a.-  -  **  II 

0 

<  afi  1 1  -  X  *  1 1  +  1 1 1  1 1  +  1 1 1  ».^-r  “’•■^  1 1  f 

<  af{  1 1 (*••■)+!  _  X  *  1 1  +  1 1  a,*->(*.<)+i  _  a;*-y(*.«)  1 1  )2 

<  af{2 1 1  -  X  *  1 1  +  1 1  x*-''(*>‘)  -  X  '  1 1  f 


a 


Then 


II  E  (B*-F'(x*)K«^^||^ 

•  €^2 

=  E  ||(B*-F'(x>,|r 

•  €^2 

<  4  E  «.•<«"  • 

«  i€Ao 


Therefore, 


I  E  (5*-F'(x*))e..e,ni^  <* 

i€/2 


This  completes  the  proof  of  (4.5). 
By  (4.4)  and  (4.5), 


||(B*-/’'(x*))(x*+‘-x‘) 

lim  - 

k^oo 


=  lim  - 

k—*oo 


||z*  +  ‘-X*|| 

I  E(^*-^'(*'))«.^V^‘-=^*)II 

•e/j 

llx*+*-x*|l 


<  lim  HE  {Bi-F'{x*))eie,^\\f  =  0  . 

.’G/ts 


')KrfMI“ 
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Theorem  4.g.  Assume  that  F,  x°  and  {x*}  satisfy  the  hypotheses  of  Theorem  4.1.  If,  for  any 
1  <  I  <  n,  «*=0  appears  consecutively  in  at  most  m  steps,  then  the  ^convergence  order  is  not 
less  than  t,  where  r  is  the  unique  positive  root  of 

=  0  .  (4.6) 

In  particular,  if  «*^0,  i  =1, ...,  n,  k  =  1,2,3,...,  then  r  =  ^  ^  1.618. 

Proof.  Notice  that  (4.3)  implies  that  there  exist  and  ^>0  such  that  1 1  1 1  ^  lor  ^ 

*  ^  ^0-  Thus,  by  Theorem  2.4, 

llx‘+^-x*|l  =  l|x‘-x*-V^’(x*)|| 

+  (|lF'(x*)-F'(x‘)lU+  ||F'(x*)-Bt||f)||x*-x*||} 

,  *-i  .  .  ^  (4.7) 

<  “  ll*‘-*'ll  +  “  Z)  -**11 

<  |o/9(  E  llx^-x*||)||x*-x*||  . 

Thus,  the  desired  result  follows  from  Ortega  and  Rheinboldt’s  Theorem  9.2.9  [8,  p.29l]. 


5.  Numerical  Results. 


We  computed  some  examples  with  tridiagonal  Jacobians  by  the  CPR  algorithm  and  Algo¬ 
rithm  2.2.  In  this  section,  we  compare  the  numerical  results  from  the  two  algorithms.  The  global 
strategy  we  used  in  computing  the  examples  is  the  line  search  with  backtracking  strategy  (see 
Dennis  and  Schnabel  [6].  If  p*  =  -F*‘*F(x*)  is  not  a  descent  direction,  then  we  try  -p*.  If  it  is 
not  a  descent  direction  either,  then  the  algorithm  fails.  The  stopping  test  we  used  is 


max 


<  e 


and  we  choose  tppx,-  = 


i<« <«  max{  I  x,*"*^'  I ,  typxi } 

10“*  and  €  =  10~®.  For  the  lower  bound  of  1  | ,  we  choose 


(5.1) 


0  =  V macheps  1 1  s  1 1  . 

We  used  double  precision,  and  the  machine  precision  is  2.22d-16. 

Example  5.1  is  new,  and  it  can  be  seen  to  be  an  extension  of  the  Rosenbrock  [9]  function 
(also  see  More',  Garbow  and  Hillstrom  [7])  to  nonlinear  system  of  equations  with  a  tridiagonal 
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structure.  Example  5.2  was  given  by  Broyden  [l]  (also  see  More',  Garbow  and  Hillstrom  [7]). 


The  results  are  shown  in  the  tables  below,  where  IT  is  the  number  of  iterations,  NF  is  the 
number  of  function  (F’(z))  evaluations,  and  LN  is  the  number  of  line  searches  in  which  the  step 
length  X<1.  ND  is  the  number  of  nondecrease  directions.  ZR  is  the  number  of  the  iterations  that 
there  exists  an  integer  j  such  that  ]  |  <  0.  *0  is  the  initial  guess. 

Example  5.1. 

/i(*)  =  8(zi-x|), 

fj{x)  =  16Xj[Xj  -  Zy_l)  -  2(1  -  Zy)  +  8(zy— Zy^l),  2,...,  H-l, 

fn{x)  =  16a'«(*»“  -  *«-i)  -  2(1  -  I.), 

n  =  9 

zl  =  (-1,  -1,  ...,  -1)^,  z2  =  (-0.5,  -0.5,  ...,  -0.5)^,  z3  =  (2,  2,  ...,  2f. 


Algorithms 

z0=zl 

IT 

NF 

LN 

ND 

ZR 

CPR 

Alg.  2.2 

22 

24 

88 

73 

15 

13 

0 

4 

0 

0 

Table  5.1  (a). 


Algorithms 

H 

II 

O 

IT 

NF 

LN 

ND 

ZR 

CPR 

AIk.  2.2 

22 

22 

88 

67 

15 

13 

0 

1 

0 

0 

Table  5.1  (b). 


Algorithms 

z0=z3 

IT 

NF 

LN 

ND 

ZR 

CPR 

AIk.  2.2 

8 

10 

32 

31 

0 

0 

0 

0 

0 

0 

Table  5.1  (c). 
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Example  5.2  (Broyden  tridiagonal  function). 
fi{x)  =  (3  -  2x,)z,  -  z,_i  -  2z,+,  +  1, 

*0  ==  *.+i  =  0, 

n  =  9, 

xl  =  (-1,  -1,  -1)^,  x2  =  (-0.3,  0.3,  ...,  -0.3,  0.3)^, 

x3  =  (-10,  -10,  ...,-10)^. 


Algorithms 

x0=xl 

IT 

NF 

LN 

ND 

ZR 

CPR 

Alg.  2.2 

5 

6 

20 

19 

0 

0 

0 

0 

0 

0 

Table  5.2  (a). 


Algorithms 

x0=x2 

IT 

NF 

LN 

ND 

ZR 

CPR 

Alg.  2.2 

5 

6 

20 

19 

0 

0 

0 

0 

0 

0 

Table  5.2  fb). 


Algorithms 

x0=x3 

IT 

NF 

LN 

ND 

ZR 

CPR 

Alg.  2.2 

8 

10 

32 

31 

0 

0 

0 

0 

0 

0 

Table  5.2  (c). 


6.  Concluding  Remarks. 

We  have  presented  an  algorithm  for  solving  sparse  nonlinear  systems  of  equations.  This 
algorithm  is  based  on  consistent  partitions  of  the  columns  of  the  Jacobians,  and  it  is  a  combina¬ 
tion  of  the  CPR-CM  algorithm  and  a  secant  algorithm.  This  algorithm  incorporates  the  advan¬ 
tages  of  the  CPR-CM  algorithm  and  secant  algorithms  in  such  a  way  that  it  reduces  by  one  the 
number  of  function  evaluations  required  by  the  CPR-CM  algorithm  at  each  iteration,  and  it  has 
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good  local  convergence  properties.  We  have  shown  that  the  SFD  algorithm  is  locally  q- 
superlinearly  convergent,  and  that  under  reasonable  assumptions,  the  r-convergence  order  of  the 

l+y/5 


SFD  algorithm  is  not  less  than 


-,  which  is  the  r-convergence  order  of  the  one  dimensional 


secant  algorithm.  Our  numerical  results  indicate  that  when  p ,  the  number  of  the  groups  in  a  par¬ 
tition  of  the  columns  of  the  Jacobian,  is  not  large,  the  SFD  algorithm  is  probably  more  efficient 
than  the  CPR-CM  algorithm. 


The  idea  exploited  here  can  also  be  used  with  Powell  and  Toint’s  [10]  work,  which  will  lead 
to  a  method  for  unconstrained  optimization  problems.  This  will  be  our  future  work. 
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